3. Example Solow model

In this jupyter notebook we will specify, solve and analyse a simple Solow model in ModelFlow.

#Required packages
import pandas as pd

# Modelflow modules
from modelclass import model
   
#for publication 
latex=0
model.widescreen()

3.1. Specify the model

We start by defining the logic of the Solow model in the Business Logic Language.

fsolow = '''\
Income          = a  * Capital**alfa * Labor **(1-alfa) 
Consumption     = (1-saving_rate)  * Income 
Investment      = Income - Consumption 
diff(Capital)   = Investment-Depreciation_rate * Capital(-1)
diff(Labor)     = Labor_growth * Labor(-1) 
Capital_intensity = Capital/Labor 
'''

3.2. Create a model class instance

After defining the Business Logic Language and storing it in the variable ‘fsolow’, we create a class instance called msolow.

msolow = model.from_eq(fsolow,modelname='Solow model')

Above the equation for Capital and Labor on the left hand side of the = (equal to) consist of an expressing diff(Capital) and diff(Labor). The equations are not normalized.

To solve a model in modelflow all equations has to be normalized. Meaning that the left hand side only consist of variables not expressions. So the function model.from_eq will normalize the equations as the first step before the model can be solved.

In this case first diff(Capital) is transformed to \(\Delta capital = capital-capital(-1)\). Then the lagged variables is moved to the right side of the =. The same goes for diff(labor).

So the normalized business language of the model now looks like:

msolow.print_model
FRML <> INCOME          = A  * CAPITAL**ALFA * LABOR **(1-ALFA)  $
FRML <> CONSUMPTION     = (1-SAVING_RATE)  * INCOME  $
FRML <> INVESTMENT      = INCOME - CONSUMPTION  $
FRML <> CAPITAL=CAPITAL(-1)+(INVESTMENT-DEPRECIATION_RATE * CAPITAL(-1))$
FRML <> LABOR=LABOR(-1)+(LABOR_GROWTH * LABOR(-1))$
FRML <> CAPITAL_INTENSITY = CAPITAL/LABOR  $

3.3. Create some data

To show what Modelflow can do, we create a Pandas dataframe with input data. And print the first 5 out of 300 observations.

N = 300  
df = pd.DataFrame({'LABOR':[100]*N,
                   'CAPITAL':[100]*N, 
                   'ALFA':[0.5]*N, 
                   'A': [1]*N, 
                   'DEPRECIATION_RATE': [0.05]*N, 
                   'LABOR_GROWTH': [0.01]*N, 
                   'SAVING_RATE':[0.05]*N})
df.head(2) #this prints out the first 5 rows of the dataframe
LABOR CAPITAL ALFA A DEPRECIATION_RATE LABOR_GROWTH SAVING_RATE
0 100 100 0.5 1 0.05 0.01 0.05
1 100 100 0.5 1 0.05 0.01 0.05

3.4. Run the model

result = msolow(df,keep='Baseline') # The model is simulated for all years possible 
result.head(29)
LABOR CAPITAL ALFA A DEPRECIATION_RATE LABOR_GROWTH SAVING_RATE CAPITAL_INTENSITY INVESTMENT CONSUMPTION INCOME
0 100.000000 100.000000 0.5 1.0 0.05 0.01 0.05 0.000000 0.000000 0.000000 0.000000
1 101.000000 100.025580 0.5 1.0 0.05 0.01 0.05 0.990352 5.025580 95.486029 100.511609
2 102.010000 100.076226 0.5 1.0 0.05 0.01 0.05 0.981043 5.051924 95.986562 101.038487
3 103.030100 100.151443 0.5 1.0 0.05 0.01 0.05 0.972060 5.079029 96.501546 101.580575
4 104.060401 100.250762 0.5 1.0 0.05 0.01 0.05 0.963390 5.106891 97.030930 102.137821
5 105.101005 100.373733 0.5 1.0 0.05 0.01 0.05 0.955022 5.135509 97.574667 102.710176
6 106.152015 100.519926 0.5 1.0 0.05 0.01 0.05 0.946943 5.164880 98.132713 103.297593
7 107.213535 100.688931 0.5 1.0 0.05 0.01 0.05 0.939144 5.195002 98.705029 103.900030
8 108.285671 100.880357 0.5 1.0 0.05 0.01 0.05 0.931613 5.225872 99.291576 104.517449
9 109.368527 101.093830 0.5 1.0 0.05 0.01 0.05 0.924341 5.257491 99.892323 105.149813
10 110.462213 101.328993 0.5 1.0 0.05 0.01 0.05 0.917318 5.289855 100.507238 105.797092
11 111.566835 101.585506 0.5 1.0 0.05 0.01 0.05 0.910535 5.322963 101.136294 106.459257
12 112.682503 101.863045 0.5 1.0 0.05 0.01 0.05 0.903983 5.356814 101.779468 107.136282
13 113.809328 102.161300 0.5 1.0 0.05 0.01 0.05 0.897653 5.391407 102.436738 107.828145
14 114.947421 102.479976 0.5 1.0 0.05 0.01 0.05 0.891538 5.426741 103.108087 108.534829
15 116.096896 102.818793 0.5 1.0 0.05 0.01 0.05 0.885629 5.462816 103.793501 109.256317
16 117.257864 103.177483 0.5 1.0 0.05 0.01 0.05 0.879920 5.499630 104.492967 109.992597
17 118.430443 103.555792 0.5 1.0 0.05 0.01 0.05 0.874402 5.537183 105.206478 110.743661
18 119.614748 103.953478 0.5 1.0 0.05 0.01 0.05 0.869069 5.575475 105.934027 111.509502
19 120.810895 104.370310 0.5 1.0 0.05 0.01 0.05 0.863915 5.614506 106.675612 112.290118
20 122.019004 104.806070 0.5 1.0 0.05 0.01 0.05 0.858932 5.654275 107.431233 113.085509
21 123.239194 105.260550 0.5 1.0 0.05 0.01 0.05 0.854116 5.694784 108.200894 113.895678
22 124.471586 105.733554 0.5 1.0 0.05 0.01 0.05 0.849459 5.736032 108.984599 114.720631
23 125.716302 106.224895 0.5 1.0 0.05 0.01 0.05 0.844957 5.778019 109.782359 115.560378
24 126.973465 106.734397 0.5 1.0 0.05 0.01 0.05 0.840604 5.820747 110.594185 116.414931
25 128.243200 107.261893 0.5 1.0 0.05 0.01 0.05 0.836394 5.864215 111.420090 117.284305
26 129.525631 107.807224 0.5 1.0 0.05 0.01 0.05 0.832323 5.908426 112.260093 118.168518
27 130.820888 108.370242 0.5 1.0 0.05 0.01 0.05 0.828386 5.953380 113.114212 119.067591
28 132.129097 108.950808 0.5 1.0 0.05 0.01 0.05 0.824578 5.999077 113.982470 119.981548

3.5. Create a scenario and run again

dfscenario = df.upd('LABOR_GROWTH + 0.002')  # create a new dataframe, increase LABOR_GROWTH by 0.002
scenario   = msolow(dfscenario,keep='Higher labor growth ') # simulate the model 

3.6. Now the results are also embedded in msolow.

  • .basedf contains the first run of the model

  • .lastdf contains the last run of the model

Also in this case the keyword keep is used. This causes the results to be stored in a dictionary msolow.keep_solutions. This can be useful when comparing several scenarios.

3.7. Inspect results

3.7.1. Using the [ ] operator

We can select the variables of interest with wildcards. This will operate the results stored in basedf and .lastdf

3.7.1.1. Look at variables starting with a C

msolow['#ENDO']

3.7.1.2. Look at all endogenous variables

msolow['labor*'].dif.plot() 
../../../_images/Example Solow_23_0.png ../../../_images/Example Solow_23_1.png

3.7.2. Using the keept solutions

As mentioned above, because the keyword keep was used. The results are also stored in a dictionary. These data can also be used for charting.

The reason for placing the results in a dictionary is to enable comparison of many scenarios, not just the first and the last.

with msolow.set_smpl(1,30):
    msolow.keep_plot('income con*' ); 
../../../_images/Example Solow_25_0.png ../../../_images/Example Solow_25_1.png
msolow.modeldash('INCOME',jupyter=1)
apprun
Dash app running on http://127.0.0.1:5001/
Dash_graph(mmodel=<
Model name                              :          Solow model 
Model structure                         :         Simultaneous 
Number of variables                     :                   11 
Number of exogeneous  variables         :                    5 
Number of endogeneous variables         :                    6 
>, pre_var='INCOME', filter=0, up=1, down=0, time_att=False, attshow=False, all=False, dashport=5001, debug=False, jupyter=1, show_trigger=False, inline=False, lag=False, threshold=0.5, growthshow=False)

3.8. More advanced topics

3.8.1. The logical stucture

Now the model has been analyzed, and the structure can be displayed.

You will find more on the logical structure here

3.8.1.1. Model structure

msolow.drawmodel( title="Model Structure", png=latex,size=(15,15))
../../../_images/Example Solow_30_0.svg

3.8.1.2. Adjacency matrix

Another way to illustrate the dependency graph is an adjacency matrix.

msolow.plotadjacency();
../../../_images/Example Solow_32_0.png

The variables [‘INVESTMENT’, ‘CONSUMPTION’, ‘CAPITAL’, ‘INCOME’] in the red area are the core of the model and has to be solved as a system.

LABOR is the prolog and can be calculated before the core is solved. While CAPITAL_INTENSE is the epilog which can be calculated after the core is solved.

Many models comform to this pattern. And for solving purpose a model is divided into a prolog, core and an epilog. Even if the core is actually consistent of several strong components.

3.8.2. The python function used to solve the model

In order to solve the model Modelflow will generate a python function which implements the model. The user will hopeful newer have to relate to the generated python code. The point of modelflow is, that the user has to relate to the specification of the business logic, not the implementation in code

print(msolow.make_los_text)
def make_los(funks=[],errorfunk=None):
    import time
    import tqdm
    from numba import jit
    from modeluserfunk import jit, recode
    from modelBLfunk import array, classfunk, clognorm, exp, gamma, inspect, jit, lifetime_credit_loss, log, logit, logit_inverse, lognorm, matrix, mv_opt, mv_opt_prop, norm, normcdf, qgamma, sqrt, sum_excel, transpose
    def prolog0(values,outvalues,row,alfa=1.0):
        try :
            pass
            values[row,0]=values[row-1,0]+(values[row,5]*values[row-1,0])
        except :
            errorfunk(values,sys.exc_info()[2].tb_lineno,overhead=9,overeq=0)
            raise
        return 
    def prolog(values,outvalues,row,alfa=1.0):
        prolog0(values,outvalues,row,alfa=alfa)
        return  
    def core0(values,outvalues,row,alfa=1.0):
        try :
            pass
            values[row,10]=values[row,3]*values[row,1]**values[row,2]*values[row,0]**(1-values[row,2])
            values[row,9]=(1-values[row,6])*values[row,10]
            values[row,8]=values[row,10]-values[row,9]
            values[row,1]=values[row-1,1]+(values[row,8]-values[row,4]*values[row-1,1])
        except :
            errorfunk(values,sys.exc_info()[2].tb_lineno,overhead=20,overeq=1)
            raise
        return 
    def core(values,outvalues,row,alfa=1.0):
        core0(values,outvalues,row,alfa=alfa)
        return  
    def epilog0(values,outvalues,row,alfa=1.0):
        try :
            pass
            values[row,7]=values[row,1]/values[row,0]
        except :
            errorfunk(values,sys.exc_info()[2].tb_lineno,overhead=34,overeq=5)
            raise
        return 
    def epilog(values,outvalues,row,alfa=1.0):
        epilog0(values,outvalues,row,alfa=alfa)
        return  
    return prolog,core,epilog